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We consider a mixture of one neutral and two oppositely charged types of molecules confined 
to a surface. Using analytical techniques and molecular dynamics simulations, we construct the 
phase diagram of the system and exhibit the coexistence between a patterned solid phase and a 
charge-dilute phase. The patterns in the solid phase arise from competition between short-range 
immiscibility and long-range electrostatic attractions between the charged species. The coexistence 
between phases leads to observations of stable patterned domains immersed in a neutral matrix 
background. 



Mixtures of cationic and anionic amphiphilic molecules 
can form thermodynamically stable structures such as 
micelles, membranes and multilamellar systems. These 
self- assembled structures have been studied as a function 
of the molar ratio of the oppositely charged molecules, 
their concentration in solution, and the ionic strength 
of the environment P, [3, H, |3l- In addition, the pres- 
ence of other neutral components leads to a large vari- 
ety of structures and the possibility of local organization 
on the surface of membranes and monolayers in emul- 
sions. Such structures are important in diverse applica- 
tions, such as the design bio-sensing devices [sj. More- 
over, they serve as model systems for the understanding 
of the properties of cell membranes where structured do- 
mains are known to be crucial to cell signaling processes 
0. In multicomponent membranes periodic nanostruc- 
tures can be immersed in a homogenous background. 
This letter discusses the possibility of this coexistence of 
periodic nanostructures with dilute phases of immicible 
cationic and anionic molecules. Besides charges on bio- 
logical membranes, the competition of short range immis- 
cibility and long range attractions leads to the formation 
of periodic structures in a large variety of systems includ- 
ing lipid mixtures JL8, 9, 10, 11], two dimensional uni- 
axial ferromagnets [l^l, reaction controled phase segre- 
gating mixtures |l3l| , and two dimensional electron gases 
in MOSFET's 

Phase separation phenomena at surfaces is well known 
[H, [H, Uzi and can be studied through simple mod- 
els of immiscibility. On the other hand, several recent 
studies have shown that mixtures of immiscible oppo- 
sitely charged molecules can form regular periodic nanos- 
tructures (or microphases) [H, [l^, [13, |2l|. To study 
the convergence of these two phenomena we use ana- 
lytic techniques as well as off-lattice molecular dynamics 
simulations of a coarse-grained model with two immis- 



cible charged molecular components. We assume that 
the membrane surface is in equilibrium and that fluc- 
tuations perpendicular to the surface are negligible. To 
begin with, we explore the phase behavior of this model 
analytically when the domains are well segregated and 
periodically ordered. We briefly describe the effect of 
different charge ratios of the molecules, but concentrate 
on the effects of the presence of the third neutral com- 
ponent and determine the conditions of coexistence be- 
tween a microstructured solid and a dilute gas of charges. 
Secondly, we describe the results of molecular dynamics 
simulations for intermediate segregation regimes. Sim- 
ulation results exhibit the expected microscopic phase 
behavior, but also indicate the limits of our theoretical 
analysis. 

We consider only two possible phases for the system 
at low temperature. One phase consists of a dense, pat- 
terned solid formed by the charged components. Its free 
energy is computed by assuming the formation of regions 
of constant particle and charge density. We ignore the 
fluctuations in the charge density and the shape of the 
interface. The second phase is homogenous, has a low 
density of charged particles, and can be treated as a two 
dimensional charged gas. Its free energy can be calcu- 
lated at low temperatures using linear response theory 
by means of the Random Phase Approximation (RPA). 
The resulting phase diagram is plotted in Fig. [T] 

Consider first the calculation of the free energy of the 
solid phase. The average absolute value of the charge 
density in this phase is it?/), (in electronic charges per 
unit area) , and the line tension between microdomains is 
7. As shown in previous work, the domains form a lamel- 
lar structure for near symmetric charge density ratios, 
while for highly asymmetric cases we have near-circular 
domains arranged in a hexagonal lattice. The periodicity 
of both of these possible structures defines a characteris- 
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FIG. 1: Phase-diagram of the oppositely charged immiscible 
mixture. Location of solid-gas coexistence curves in p — T 
plane is shown for several different values of the Bjerrum 
length, Ib- Increasing the Bjerrum length requires higher val- 
ues of short range interactions e for phase coexistence. The 
curves shown for the left branch of the coexistence region lie 
in the shaded region of the sketched full diagram in the cen- 
tral inset. The solid phase has constant density psoUd- Fig- 
ure refers to a symmetrically charged mixture that generates 
lamellar structures in the solid phase. 



tic length L. Then, the free energy contribution per area 
A, due to the microphase formation has the form 



Here, the Bjerrum length is Ib 



(1) 



where fc^ is 



the Boltzmann constant and T is the temperature, and 
£ is the dielectric permittivity of the medium. The co- 
efficients si and S2 are dimensionless quantities that de- 
pend on specific shape of the microdomains. si is the 
ratio of the microdomain interface length within a unit 
cell to the size of the cell. S2 is the integral of the di- 
mensionless Coulombic potential over the whole space, 
averaged over a unit cell. These values are calculated 
explicitly in Ref.[l^ for both the lamellar and hexago- 
nal cases. Minimizing the free energy with respect to 
the size L, we obtain an optimal characteristic length 
Lrmn = {si/ 82)^^^ Lq, with Lq = (7/('s ) )^/^ • Evalua- 
tion of the free energy density at that point results on a 
value of 2(siS2)i/2/o/i§, where /o = hh^Y^^ MM- 
The net free energy of the solid phase also includes the co- 
hesive energy that arises from segregation of the charged 
molecules from the neutral component. To be able to 
compare energies with the more dilute ionic gas, and with 
results of the molecular dynamics simulations, it is nec- 
essary to introduce a model of the individual molecules. 
In the simulation, their Wan der Waals interactions are 



r < Tr 
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where the potential is cut at a radius Tc = 2.50' for 
oppositely charged molecules and = 26 for similarly 
charged molecules, where a is an effective molecular ra- 
dius. This selection of cut-offs produces a net immis- 
cibility between species of magnitude e per contact be- 
tween oppositely charged particles. In our description of 
the solid phase, at low temperatures, we simply assume 
a constant density given by the hep packing of spherical 
molecules of radius a. For this close packing arrangement 
we have the constant density PsoUd = ■ The effec- 
tive cohesive energy per unit area can then be written 
in terms of e, with 6 favorable contacts with neighbors. 
Using a similar approximation, the effective line tension 
is ^/6e/a. Inclusion of the cohesive energy, leads the final 
result for the solid phase as: 
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Next, consider the free energy of the gas phase. RPA 
[2^ is used to calculate the free energy as a function of 
the relative strength of the short range attraction and 
the electrostatic interactions. The method involves an 
expansion of the free energy of the system in terms of 
density fluctuations, neglecting all terms of larger than 
second order. For a general system of interacting i,j 
charged monomers, the partition function is written as 
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Here, pi represents the density of the ith component and 
Pk represents the Fourier transform of the component 
densities. Zq includes the k = Q and self energy terms. 
U^-' is the sum of the interaction energies, consisting of 
the short range interactions due to the excluded volume 
and hydrophobic interaction, as well as the long range 
contributions due to the electrostatic energy. The elec- 
trostatic contribution to the internal energy matrix uses 
the two dimensional Fourier transform of the screened 
Coulomb interaction, iJ^] = 2TrziZjlB{k^ + k'^)^^^^. The 
free energy is then 
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The first sum (the fc = terms) vanishes due to charge 
neutrality. The second term is a standard entropic term. 
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The third term is the electrostatic contribution due to 
density fluctuations, which is calculated by integrating 
over the possible values of k from to 27r/a, where a 
is the molecular size. We neglect terms arising from ex- 
cluded volume interactions, as the calculation is domi- 
nated by the electrostatic terms. We only consider the 
limit of no screening due to ions in the surrounding solu- 
tion (matching the conditions of the calculation for the 
solid phase and for the molecular dynamics simulations). 
The electrostatic contribution is found to be, 
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where kf^ — 2'kIb{p+z'\. + p+zj_). The above equa- 
tions simplify when considering the charge neutrality 
constraint, 2+p+ = Z-P-. The total free energy per unit 
area for the gas phase, in terms of p where p = p+ + p-, 
is then 
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potential from < r < cr and a classic 6-12 Lennard 
Jones potential from cr < r < 2.5(t. 

To determine the phase coexistence of the solid, pat- 
terned phase and the dilute gas, we use the common tan- 
gent rule. Since the density of the solid phase is assumed 
fixed, the equation 
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is solved for the gas phase density p. We plot the two 
phase coexistence line calculated in this way, for low tem- 
peratures, in Fig. 1. The figure refers to the case of equal 
charge density — 1. The phase diagram shape de- 

pends of course on the relative strengths of the Coulomb 
interaction (through the charge density) and on cohe- 
sive/immiscibility parameter e. With increasing strength 
of the electrostatics, the transition points occurs at lower 
temperatures, while increasing values of the short range 
interaction increase those temperatures. At higher val- 
ues of the density, nonlinear corrections to RPA including 
short range correlations and ion association 0, [1^. Fig. 
1 sketches an extrapolation of the results to higher tem- 
peratures when the solid phase retains its near constant 
density. 

To explore higher temperatures and the effects of fluc- 
tuations, we carried out a set of molecular dynamics 
simulations. The model systems used are composed, in 
the symmetric case, of a mix of = 1000 positively 




FIG. 2: Simulation snapshots for systems with charge ratios, 
a-b) z+/z- = 1/3 and c-f) z+/z- = 1. Frames a-b illus- 
trate hexagonal order for different densities pa = 0.10, pi, = 
0.36. Frames c-d show the effect of increased immiscibil- 
ity Ec = SAksT, ea = SAksT, for a fixed Bjerrum length, 
Ib/(J = 0.5. Frames e-f show the transition from a homoge- 
neous microphase to a phase segregated state for te = 2.0kBT, 
6/ = 2.SkBT. 



charged and = 1000 negatively charged units in a 
simulation box of size D^ with D = 66ct. In the asym- 
metric case we used — 900 and N^^ = 300 units 
with charge (-3). The molecules are confined to a plane 
perpendicular to the Z axis, with periodic boundary con- 
ditions in the X and Y directions. The molecular dy- 
namics simulations, with constant N,V,T were performed 
using the Espresso simulation code of the MPIP-Mainz 
group of Polymer Theory and Simulation [2^. We ex- 
plored regions of the phase diagram at surface densities of 

P = ''^^^4%2 ^^'^ = 0.36. The potential between charges 
is a Coulomb potential, Uc — ^s2l3ill2.^ Further param- 



eters are as described in detail in previously published 
work p^ . Late-time snapshots of the simulations are 
shown in Fig. 2. To compare the behavior observed in 
the simulations to our theoretical results, we point out 
that the simplest extrapolation of the results from the 
strong segregation limit to higher temperatures, consist 
on identifying the characteristic length scale of system 
(the size of microdomains) , up to a constant, with the 
form L = (7/(/biA'))'/' ~ (e/(^B^3/2))i/2_ use this 
expression for extrapolation, we abandon the assumption 
of a solid phase density fixed by the molecular radius, but 
take ip as the density determined by the simulation at fi- 
nite temperatures. All our simulations exhibit the same 
expected basic behavior. At small values of e or low Bjer- 
rum lengths (high temperature), positive and negative 
regions develop on the surface and as the temperature 
decreases the domains increase in size. In all the sim- 
ulations, the individual molecular components exhibit a 
hexagonally close-packed structure (as was assumed in 
our analytical approach), with density fiuctuations de- 
pendent on the temperature. For asymmetric charge ra- 
tios, — = 1/3, we show in Fig. 2a the formation of a 
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FIG. 3: The location of the peak k* in the structure factor 
S{k) as a function of e/(Zs?/'^'^^). The hnear fit shows agree- 
ment with the scahng predicted by strong segregation theory 
of power —0.5. 

hexagonally patterned 'island' at p = 0.10. For larger 
densities, as in Fig. 2b the solid phase occupies a larger 
fraction of the space, but exhibits more clearly the order- 
ing. 

Next, as shown in Fig. 2c-d, for symmetric charge ra- 
tios, -J-*- = 1 the preferred microstructure is lamellar. At 
the values of the parameters used, we observe as well the 
phase separation between solid patterned and neutral re- 
gions. The temperature decrease from snapshot c to d, 
clearly modifies the fluctuations of both types of inter- 
faces: between charged regions, and between the solid 
and neutral phases. At lower temperatures (Fig. 2d) 
the interfaces are much sharper and exhibit lower shape 
fluctuations. In these simulations, a note and an inter- 
esting feature of the solid-gas interface: the orientation 
of the lamella is perpendicular to the interface. While 
the charged domains have symmetric interactions with 
the neutral region, the alignment must be a result of 
minimization of the local electrostatic energy. To some 
extent, this feature also appears in the asymmetric case. 
Fig. 2a-b, where both charged domains appear the at the 
interface wit the neutral domain. 

For symmetric charge ratios, but for weaker electro- 



static interactions. Fig. 2e-f show the transition from a 
solid to a solid-gas coexistence phase. At low values of the 
cohesive energy e the system shows the lamellar pattern- 
ing but has large voids between the charged domains; the 
neutral component is attracted to the interfaces where it 
can reduce the effective line tension between domains. 
On further increase in the cohesive parameter, the coex- 
istence region is reached, and the neutral regions segre- 
gate to form their own phase, as shown in Fig. 2f. The 
lower values of the Bjerrum length in these cases produce 
larger lamellar sizes, compared with those of Fig. 2c-d. 

To study the behavior of the domain size (lamellar 
width) in the symmetric case, we have determined the 
structure factor S{k) = (pkP-k), the Fourier transform 
of the correlation function of the charged components, 
for the late stages of a number of simulations with differ- 
ent values of the interaction parameters. S{k) displays 
a peak at values k* corresponding to the inverse lamel- 
lar spacing in the direction perpendicular to the lamellas, 
and thus must scale as k* ~ (e/(ZB?/''^^^))~^^^- Indeed, in 
our simulations, we find that the position of the structure 
factor peak can be fit through a line of slope —.47-1-/ — .02 
when plotted against the group e/(lB'0^^^)i ^is shown in 
Fig. 3. 

By use of analytic techniques, combining RPA and 
strong segregation theory, as well as off-lattice molecu- 
lar dynamics simulations we have demonstrated the clear 
possibility of coexistence of structured regions, made of 
charged components, in a matrix rich in neutral compo- 
nents. The competition between immiscibility tendencies 
between molecules, against the attractive interactions be- 
tween oppositely charged species, provides a method to 
generate well-controlled, self-assembled surface patterns. 
Furthermore, suitable extensions of the results from the 
strong segregation limit to higher temperatures provide 
a clear way to quantitatively address the properties of 
these materials in practical regimes. 
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